
SEPARATED FLOWS NEAR THE NOSE OF A 


BODY OF REVOLUTION 


Supported by NASA Grant NCC2-280 

{NASA-CR- 17710?) SEPARATED FLOSS NEAR THE N86-30091 

NOSE OF A BODY OF REVOLUTION Final Project 

Report, 1 Jan. ,1964 - 28 Feb. ,1986 {Clarkson 

Univ.), 89 p CSCL 20D Unclas 

G3/34 43217 


Submitted by 


Dr. 


S.P. 


Lin, 


Principal 


Investigator 


Department of Mechanical and Industrial Engineering 


r 


Clarkson University 


Potsdam, NY 13676 


c 


TABLE OF CONTENTS 


PART I - Project Identification Information 
PART II — Summary of Completed Project 
PART III - Technical Information 


APPENDIX 



PART I - PROJECT IDENTIFICATION INFORMATION 

1. Principal Investigator: Dr. S.P. Lin, Professor, Mechanical and Industrial 

Engineering 

2. Institution: Clarkson University, Division of Research, Potsdam, NY 13676 

3. NASA Program: NASA-Ames Cooperative Program, NCC2-280 

4. Award Period: From 1/1/84 to 2/28/86 

5. Cumulative Award Amount: $80,063 

6. Project Title: Separated Flows Near the Nose of a Body of Revolution 

7. Technical Program Officer: Dr. Gary Chapman, NASA-Ames 


1 



PART II - SUMMARY OF COMPLETED PROJECT 

The solution of the Navier-Stokes equations for the problem of cross-flow 
separation about a deforming cylinder is achieved by a novel method of 
iteration. It is shown that the separation starts at the rear stagnation point 
and the point of primary separation moves upstream along the cylinder surface. 
The solution clearly indicates how the secondary separation may follow in time 
the primary separation. By invoking the analogy between the cross-flow 
structure around a three-dimensional body of revolution and the two-dimensional 
flow around a deforming cylinder, we postulate that the flow separation does not 
originate at the nose-tip. We also postulate how the open separation may 
originate near the nose-tip. 

A general method of linear stability analysis for non-parallel external 
flows has been constructed. It consists of representing the eigenfunctions with 
complete orthogonal sets and forming characteristic equations with the Galerkin 
method. The method is applied to the Kovasznay flow which is an exact solution 
of the Navier-Stokes equation. The results show that when the critical 
parameter is exceeded, there are only a few isolated unstable 
eigen-frequencies. Hence the concept of neutral curve which is useful for 
parallel flows becomes irrelevant. Another exact solution, i.e. Taylor's vortex 
array is shown to be absolutely and monotonically stable with respect to 
infinitesimal disturbances of all frequencies. The flow is also globally, 
asymptotically and monotonically stable in the mean with respect to 
three-dimensional disturbances. This result forms the sound foundation of 
rigorous stability analysis for non-parallel flows, and provides an invaluable 
test ground for future studies of non-parallel flows in which the basic states 
do not possess exact solutions. The application of this method to the study of 
the formation of spiral vorticies near the nose of a rotating body of revolution 
is underway. The same method will be applied to the stability analysis of the 


J 
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reversed flow over a plate with suction. This flow is an exact solution of the 
Navier-Stokes equations, and was obtained by us. 
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PART III - TECHNICAL INFORMATION 


1. Scientific Collaborators 

The following individuals have collaborated with the principal 
investigator. Dr. S.P. Lin on this project: 

David Mekala - Graduate Research Assistant, completed his M. S. thesis in 
June 1985 

H.B. Chen - Graduate Research Assistant, completed his M. S. thesis in 
January 1986 

K. El-Rais - Ph.D. Candidate 

E. Ibrahim - Ph.D. Graduate Research Assistant 

Prof. R. Hessel and Prof. D. Valentine of Clarkson University have been 
collaborating with the principal investigator on some numerical aspects of 
the proj ect 

M. Tobak and G. Chapman of NASA-Ames have provided valuable physical 
insight and scrutinized some technical detail. 

2. Abstracts of Theses 

A. "Migration of the Separation Point on a Deforming Cylinder" 

M. S. Thesis by David Mekala 

An iterative scheme of solving the Navier-Stokes equations for the 
two-dimensional cross-flow about a deforming cylinder is given. In the 
first approximation, the effect of vorticity convection is neglected. 
This effect is taken into account by successive iterations. Thus in 
the k-th iteration, the convected vorticity field appears as the source 
term in an equation of transient diffusion of vorticity. A novel 
integral transform is used to reduce this transient vorticity equation 
into a k-dimensional heat equation. The bounded solution of this 
equation with the appropriate boundary conditions is given. Based on 
analytically obtained flow field, the equation which gives the 
locations of zero shear stress at the surface of the deforming cylinder 
is found. The final solution for the separation angle is expressed as 
the sum of integrals involving, explicitly, the Reynolds number, and, 
implicitly the dimensionless cylinder deformation. Gauss-quadrature is 
used to evaluate these integrals for numerical evaluation. The 
specific forms of the cylinder motion used in the present numerical 
calculations include impulsively started uniform motion about a rigid 
cylinder, motion consisting of constant acceleration about a deforming 
cylinder and motion consisting of impulsively started uniform motion 
about deforming and contracting cylinders. The general solution is 
applied at the surface of the cylinder to determine the locations of 
primary separation points based on the first two iterations. 
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B. "Stability of a Non-Parallel Flow," M. S. Thesis by H.B. Chen 

There have been many brilliant works on the stability of parallel flows 
during the past decades. Unlike the stability analysis of parallel 
flows, the stability analysis of non-parallel flows is not yet firmly 
established. Among other difficulties, the stability analysis of 
non-parallel flows suffers from some serious uncertainty associated 
with the approximate nature of the basic state of the non-parallel 
flows. Kovasznay obtained an exact solution to the Navier-Stokes 
equations for a non-parallel flow. He suggested that this flow may be 
used to describe the flow downstream of a two dimensional grid. We 
take this full non-parallel flow as our basic flow of stability 
analysis. The linear stability of the Kovasznay flow is analyzed. The 
solution to the obtained governing equation of the linear stability is 
expanded in two complete sets of orthogonal functions which satisfy the 
appropriate boundary conditions. The characteristic equation is then 
obtained by the Galerkin method. The stability of the flow is 
characterized by the Reynolds number R and another parameter B 
representing the length of the closed wake. It is shown that the 
Kovasznay flow is absolutely stable when B = 0 for any R. The 
numerical results show that there exists a critical B below which the 
flow is stable for a given R. Above the critical B, the flow becomes 
unstable only with respect to finite number of discrete frequencies of 
infinitesimal distrubances. 

3. Publication Citations 

a. S.P. Lin and M. Tobak, "Reversed Flow Above a Plate with Suction," AIAA 
Journal 24^ 334, 1986. . 

b. S.P. Lin, D. Mekala, G. Chapman and M. Tobak, "Migration of the 
Separation Point on a Deforming Cylinder," submitted to AIAA Journal. 

c. S.P. Lin and M. Tobak, "Spectral Stability of Taylor's Vortex Array," 
to appear in Physics of Fluids. 

d. S.P. Lin and M. Tobak, "Stability of Taylor's Vortex Array," submitted 
to J. of Fluid Mech. 

e. S.P. Lin, H.B. Chen and M. Tobak, "Stability of Non-Parallel Flows," 
submitted to Physical Review Letters. 

The technical details of the above works are given in the Appendix. 

4. Continuing Work 

Building upon the results of the completed work listed above, we will 

complete the following unfinished works. 

a. "Separated Flows about a rotating Body" (see Appendix for detailed 

description) 

b. "Formation of Secondary Vorticies on a Circular Cylinder." 
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APPENDIX 


TECHNICAL DESCRIPTION 

* Instability of Non-Parallel Flows 

* Spectral Theory of Taylor* s Vortex Array 

* Stability of Taylor* s Vortex Array 

* Reversed Flow Above a Plate with Suction 

* Separated Flow about a Rotating Body 

* Migration of the Separation Point on a Deforming Cylinder 
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Instability of Non-parallel Flows 

S. P. Lin and H. B. Chen 
Clarkson University, Potsdam, NY 13676 

and 

M. Tobak 

NASA Ames Research Center 

Moffett Field, CA 94035 

This letter gives a general method of stability analysis of non-parallel 
flows which varies in the direction of flows. The solution to the governing 
partial differential equation is constructed with two complete orthogonal 
sets. The characteristic equation of disturbances is then obtained with the 
Galerkin method. The modified Kovasznay non-parallel flow which is an exact 
solution of the Navier-Stokes equation is used to demonstrate the method. 1 The 
results show that there are countably infinite numbers of stable eigenmoaes 
but only a few isolated unstable modes. Unlike the parallel flows, a neutral 
curve does not exist. 
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Stability analysis for parallel flows is quite well developed, as can be 
seen from recent books on the subject^ A far less developed analysis is 

that for non-parallel flows which are encountered in almost every external 
flow. The importance of the non-parallel flow effect even on almost parallel 
flows are increasingly recognized^ In linear stability analysis of 

parallel flows, one obtains the neutral curve which divides the parameter 
space into stable and unstable regions. In such an analysis, the disturbance 
may be represented by the Fourier integral, and the governing equation of 
stability is reduced to an ordinary differential equation. This implies that 
the wave length of the Fourier component of disturbances may be varied 
continuously. Consequently, a continuous neutral curve may be obtained in the 
wave number vs. flow parameter space. In weakly non-parallel flows, the 
concept of a neutral curve is still applicable. However, in truely 
non-parallel flows the basic flow varies in the flow direction as well as in 
other directions, and the disturbance cannot be decomposed into independent 
Fourier components with arbitrarily close neighboring wave lengths. As a 
consequence, a continuous neutral curve does not exist in general. In fact, 
the disturbance must satisfy the partial differential equation which may admit 
only isolated eigen-modes. A general method of stability analysis for 
non-parallel flows has not yet appeared. Existing stability analysis for 
weakly non-parallel flows appear to be all applied to the basic states which 
can be described only approximately. While these studies are motivated by 
practical needs, they entail the theoretical disadvantage of having to tailor 
the method of analysis to fit the degree of approximation involved in the 
basic state. Thus, if we are t^o seek a general method of linear stability 
analysis for non-parallel flows, it may be more advantageous to try it on a 
non-parallel basic state which is an exact solution of the Navier-Stokes 
equation. This will free us from concern, over the possible inadequency of the 
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approximation involved in the basic state. Motivated by this thought, we 

apply our general method of linear stability analysis to the modified 
12 

Kovasznay wake flow which is an exact solution of the Navier-Stokes 

equation. Our method consists of constructing the disturbance with two 

complete orthogonal sets and solving the characteristic equation obtained by 

13 

use of the Galerkin method. A similar method has been given by Orszag for 
internal flows . 

Consider two-dimensional flows of an incompressible Newtonian fluid. The 
X and Z Cartesian components of the velocity field are related respectively to 
the Stokes stream function Y by 

u -v v = -v 

where subscripts X and Z denote partial differentiations. In terms of Y, the 
governing equation is 

(3 t -v7 2 )V 2 T = '? x V 2 ¥ z - T z 7 2 t x , (1) 

2 

where t is time, v the kinematic viscosity# V is the Laplacian operator, the 
subscripts X and Z denote partial differentiations with respect to the space 
variables. By use of the following dimensionless variables 
**= Y/(Ud) , ( x , z ) = <X,Z)/(d) , 

V— V(d) , r = t/(d/U) , 

where U and d are respectively the characteristic velocity and length, (1) can 
be rewritten as 

O - I 7 2 )v 2 * = * V 2 * - * V 2 * . (2) 

l IV X Z Z X 

whore R = Dd/v is the Reynolds number. 

Kovasznay found the following exact solution of (2) 

= z - B sin (2uz) exp(Mx) , (3) 
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where 


M = 3*( R-/ R 2 +16tt 2 ) < 0 


He chose B to be 1/2 tt so that the stagnation points are located at x=o 

and z=+n, n being any integers. Here, we chose B to be a free parameter, and 

consider only the region x>o so that the velocity is bounded everywhere. The 

stream line pattern of this flow for R=4 0 and B=lo.O is given in Fig. 1. 

Note that the length of the closed wake, i.e. the distance between r=0 and the 
stagnation point is given by 


L M 2ttB ’ 

Since M<0, for the existence of the closed wake we must have 2 jtB> 1. Thus, B 
is an independent parameter which represents the strength of suction and blow 
at x=0. When B=0, the flow is a uniform stream of speed U in the 
x-direction. When 0<BO./2n, the flow is sucked around z=n toward x=0 but 
blown away from x=0 around z=2n-l ; however it is not large enough to form a 
closed wake. Kovasznay suggested that this solution can be used to describe 
the wake flow behind a grid of wires with uniform spacing d. 

The stability of this modified Kovasznay flow is now analysed. Let the 
perturbed stream function $ be written as 


* = * + 4> , (4) 

where <p is the stream function perturbation. Substituting (4) into (2) and 
neglecting nonlinear terms, we have 

O -|-V 2 )7 2 $ = 4 7 2 * - 4 V 2 <p + Mf(« -1)4 - * 4 ] (5) 

T K X Z Z X z x xz 

It follows from (3) that 

( 

^ ** ^ B sin (2 ttz) and ij; = 1 -$> B cos (2wz) s (6) 

x xo z zo 

where 

4 xq = -MF(x), 4 zq = 2wF(x) , F(x) = exp (Mx) 
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The coefficients in (5) are functions of x and z. Thus, it cannot be 
reduced to an ordinary differential equation by Fourier decomposition as is 
customarily done in the parallel flow stability analysis. We seek the 

solution of (5) in the form 

N N 

♦ = £ Cn(x,x) sin 2nuz + £ Dn(x ,x) cos 2mrz 

n=l n=l 

The solution is expected to approach the exact solution as N***. The chosen 
z-dependence in the above solution ensures that the stability behavior is the 
same in each wake. Physically, this implies that the disturbances are 
distributed randomly without preference of any particular wake. Substituting 
this assumed form of solution into (5), we have 
12 2 

(3 - — E )E (Cn sin 2mrz + Dn cos 2mrz) 
t R n n 

— 2 

= 2mnj> B sin 2 tt z (E - M)(Cn cos 2mrz - Dn sin 2mrz) 
xo n 

— • 2 

- [M $ B cos 2irz + (1 - ^ cos 2 ttz)E ](Cn sin 2nirz 
zo zo n 


+ Dn cos 2mrz) , 
x 


(7) 


where 


E 2 = d - (2nn) 2 . 


n xx 

The solution of (7) will be achieved by use of the Galerkin method. 
Multiplying (7) by each member of the orthogonal set sin 2mnz, m=l to N, then 
integrating over one period in z and demanding the weighted residual over this 
period to be zero, we have 


“Vi 


E 2 D]C = + B[w(m-lH (E 2 - -M) + % 4 > (E 2 ,tM)D]C . 

m m xo m-1 zo m-1 l m-1 




( 8 ) 


where D = 3 and N > m > 1. Note C =0. 
x ~ ~ o 
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Similarly, multiplying (10) by cos 2mnz and integrating over one period 


in z, we have 


[O -A E z ) E + E D]D = D - Tr(nrfl)ii> ](E - M)D 

t R m m mm zo xo m+1 


nrfl 


+ Blh <J> D- Tr(m-l) <|> ] (E 2 -M)D J (9) 

zo xo m-1 m-1 

where N _> m >0, and D. = 0 if j <0. 

Thus the odd mode and even mode disturbances are decoupled. The odd mode 

solution can be obtained from (8) which is a set of N linear fourth order 

ordinary differential equations in N unknowns C (m=l to N) . (9) is a set of 

m 

(N+l) linear fourth order ordinary differential equations for the even mode 
solution (m=0 to N) . Note that the nr-th Fourier mode is weakly coupled 

with its two immediate neighboring modes in (8) and (9). 

The instability of the basic flow with respect to the even mode 
disturbance will be investigated first. We assume that the ultimate stability 
or instability does not depend on the initial transience but depends on the 
property of the discrete modes 


D 

m 


K o lr 

E exp(aiT)^. (x) B 
k=0 k m 


( 10 ) 


« k s 

where to is the complex eigenr-frequency,B ib-.(x) is the complex disturbance 

mk 

amplitude which varies with x. We assume that the disturbance does not occur 
at x=0 but vanishes at infinity, and thus chose the representation 


^ = [x S expf-xJJ^L^CxJ/fkirCs+k+l)]^ , 

8 14 

where I* is the Gama function, and L^x) is the modified Lagurre polynomial 

s 

of order k. $ is chosen to be 4 in this analysis. Note that ^ form a 
complete orthonormal set of functions which satisfy the vanishing boundary 
conditions at x=*0 and x-h». The eigen-frequency and the amplitude coefficient 

k 

B^will be determined by use of the Galerkin method. Substituting (10) into 
(8), we have 
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where the summation sign over k is omitted, and 

L . * (q»R-E 2 )E 2 + RE 2 D 

ml mm m 

L = ir(m-l)^ (E 2 ,-M) + hty (E 2 -M)D 
m2 xo m-1 zo m-1 

L o ~ (E 2 -M)D - ir(cri-l)* (E 2 , -M) 

mJ zo nr+1 xo nrH 

Multiplying (11) with (1=0 to K), and integrating over 0<x< we have 


Q i + R i.» B + s ,.. B * . 

mkl m-1 mkl m mkl urfl 


- RwR , .B = 0 
mkl m 


( 12 ) 


where 


Vi " BR I *1 L J*l Ux 


R - ■ 1 0*1 L ml ( *k )d, ‘ 


mkl 


R 


mkl 


■ t 0 *K ( K )dx 


(13) 


S mkl ' BR ^ • 


o r l m3' r k' 


It can be seen from (12) that the Fourier modes (m) are weakly coupled 
but the Lagurre modes (k) are strongly coupled. The system (12) can be 
written in a matrix form as 


{ CG3 + R » i CH]}(V) = 0 


(14) 


where 


[G] - 


R , a S 
okl ok! 


^mki R mkl S mkl 


< Nkl 


*1^1 


[H] = 


R 


okl 


R 


mkl 


R Hkl 


and (V) is the eigen-vector defined by 

,„.T ,_k _k „k k k _k _k» 

(V) (B 0 » B i "• B m -l ’ B m * B n,+1 •" B H-1* B N 5 


in which the superscript T denotes transpose, and in each m-th mode k ranges 
from 0 to K. (14) is a system of (K+l)xN equations in (K+l)x(N) unknown 



k lr 

eigen-vector components- Note B “ 0, and are truncated to form a 

determinate system (14). The non-trivial solution of (14) exists only if the 

determinant of its coefficient matrix vanishes, i.e. 

| [G] + w.R[H]| =0 . (15) 

The solution of (15) was achieved by use of the IMSL routine EIGZF. The 

required integrals given in (13) were evaluated with the adaptive Gauss 

quadrature^ . All computations have been carried out with double precision on 

the Gould PN9780. Table I gives the complex eigen-frequencies of unstable 

modes for which the real parts are positive. The rest of eigen-frequencies, 

which are not listed, have negative real parts. At R=1000, no unstable modes 

with finite frequencies exist for values of B smaller than 6. Hence for 

R=1000, the critical value of B below which the flow remains stable with 

respect to the even mode of disturbances is 6. However, referring to Table I, 

we note that when the flow becomes unstable it does so only at a few isolated 

frequencies. This is quite different from the situation in parallel flows. 

When the flow parameter exceeds the critical value in a parallel flow, the 

flow becomes unstable with respect to disturbances of infinitely many 

frequencies within a given band-width. Comparing the eigenvalues in Table 1 

} 

obtained with two different sets of K and N, we see that K=N=10 is not yet 
large enough to yield accurate eigenvalues. However, the convergence of the 
method is expected. The main point in this letter is to point the irrelevance 
of the notion of neutral curve in non-parallel flows, and that there are only 
a few eigen frequencies at the onset of instability. 

The stability analysis for odd mode can be carried out along the same 
line starting with (8). It turns out that the flow is more stable with 
respect to the odd mode disturbances. Exact eigenvalues for both modes will 

be obtained with more extensive computation and reported elsewhere in the near 
future . 

This work was supported in part by a NASA Grant CG-280. 
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TABLE I. Unstable modes, R=1000. 


B 

N.K 

h) 

0)* 



r 

1 


9 

. 

.• — 

6 

10 

0.16 

10.4 


9 

0.10 

10.0 

8 

10 

0.69 

13.8 


9 

0.60 

12.3 

10 

10 

1.18 

17.3 
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Spectral stability of Taylor's vortex array 

S, P. Lm 

Mechanical and Industrial Engineering Department 
Clarkson University, Potsdam, NY 13676 
M. Tobak 

NASA Ames Research Center 
Moffett Field, CA 94035 

Abstract 

It is proved that the two-dimensional Taylor vortex array, which is an 
exact solution of the Navier-Stokes equation, is absolutely and monotoni cal ly 
stable with respect to infinitesimal disturbances of all discrete frequencies. 

PACS (47,51,44) 
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Consider the two-dimensional flows of an incompressible Newtonian fluid 


The governing Navier-Stokes equations written in terms of the Stokes stream 
function is 


7} eta 


\ 2 / 

= 


(i) 


* mw "A * 

nu 

where t is time, (X, Z) are the Cartesian coordinates, v is the kinematic 
2 

viscosity, ^ is the Laplacian operator, the subscripts stand for partial 
differentiation, and is the sole component of the vorticity perpendicular to 
the X-Z plane of flow, i.e. 


= . 


( 2 ) 


X 

T(k‘A. 


In terms of the following dimensionless variables 


t » t/(d 2 /v) , (X,Z) » (x,z)/d , ip » f/ v , 5 = V(v/<i 2 ) 


where d is the characteristic length, ( 1 ) can be written as 


{ pfc: Uf c ) 

i -pu 


? + Ip £ - & Z = V r . 

x z x y 7 l z s 


(3) 


S - V . 


Taylor obtained an exact solution of (3) given by 


ip = A cos itx cos nz exp(-2n x) . 


(4) 


rr 


where A is a constant representing the strength of the flow. This stream 



function represents a time dependent two-dimensional array of counter-rotating 
vorticies enclosed in squares of a dimensional side-length d. Note that 
r = (K+l)/2 and z = (K+l)/2 are stream lines, where K is any integers. Note 
also that the x- and z-components of the velocity are respectively given by 


- ip - -Ancosnx sinnz exp(-2n x) 
z 


and 


v - Ajtsinrrx cosjtz exp(-2n t) 


Hence the above mentioned two sets of stream lines intersect at saddle points, 
and the centers of vorticies are located at the stagnation points (x,z) = (M,N) 
where M and N are any integers. 

The stability analysis of this flow does not seem to have yet appeared. 

Most of the exact solutions of the Navier-Stokes equations represent parallel 

flows. The linear stability of parallel flows is governed by the well-known Orr 

2 3 

- Somme rf eld type ordinary differential equation. ' The linear stability of the 
non-parallel flow described by (4) is governed by the partial differential 
equation 





*;«z 


= 7 


(5) 


where £= V ip , and ip ' and are respectively the perturbation stream function 
and vorticity. (5) is obtained by substituting 

ip = ip + ip ’ and £ = z, + £ ' 


3 



into (3) and retaining only the linear terms, 
chara c teri s t i c frequency w, i.e. 


Con si de r/. di s t urbance s with the 


u3 

<f> ' = ♦e“ t ; C ' = £e WT , { = ? 2 <?> . 

V 




It is also assumed that the perturbation vanishes at infinity. Substitution of 
(6) into (5) yields 




(<f> §) - ((pc) - V 5 = 0 

z x x z 


(7) 


Multiplying (7) with £♦, where the upper star denotes the 
complex conjugate, and then taking the sum of the resulting equation with its 
complex conjugate, we have 


2w.5S*+ (tV(*t) v - 5> + frUc) - > (8) 

/iy z x xz zx xz 


+ {c. c } = 0 , 


where cu 

r 


is the real part of w. Integration of (8) gives 


< 2 to = <-<* (§ 5 * 

r z 


+ 5 5 *) + * ( 54 * 

X X X z 




- 5<{* <> 
* X z 


+ U) + 

X z 




.) * V. V£ - £7 .7?* 


(9) 


where < .> signifies integration over the infinite domain. In arriving at (9), 
integration by parts was applied to every term in (8) except the first one, and 
then the condition that the perturbation vanishes at infinity was invoked. 
Integrating by parts once again, and invoking the divergence theorem we have from (9) 


4 


where ^ is the absolute value of Hence 


0) 

r 


<lvd 2 > 

<UI 2 > 


(id 


This completes the proof that Taylor's vortex array is absolutely and 
monotonically stable with respect to infinitesimal disturbances of all 
frequencies as long as the viscosity is positive. The last condition follows 
from the definition of T and (6). 

The convective diffusion of perturbed vorticity happens to cancel out in 
Taylor's vortex array to yield (11). The perturbed flow is dominantly 
dissipative. However, this is not obvious a priori. It should be pointed out 
that the stability we have proved is the large time asymptotic stability of the 
discrete spectrum. The extension of the notion of asymptotic stability of 
steady basic flow to the unsteady mean flow demonstrated in this note does not 
involve any physical or mathematical ambiguity. 

ACKNOWLEDGMENTS 
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Stability of Taylor's Vortex Array 
S.P. Lin 


Mechanical and Industrial Engineering Department 
Clarkson University, Potsdam, NY 13676 
M* Tobak 

NASA Ames Research Center 
Moffett Field, CA 94035 

Abstract 

It is proved that the two-dimensional Taylor vortex array, which is an 
exact solution of the Navier-S tokes equation, is globally, asymtot ically 
and monotonically stable in the mean with respect to three-dimensional 


disturbances . 



2 


Consider the two-dimensional flows of an incompressible Newtonian fluid. 
The governing Navier-Stokes equations written in terms of the Stokes stream 
function is 


"t + Vx - Vz • ^ 


( 1 ) 


where t is time, (X, Z) are the Cartesian coordinates, v is the kinematic 
2 

viscosity, V is the Laplacian operator, the subscripts stand for partial 
differentiation, and n is the sole component of the vorticity perpendicular 
to the X-Z plane of flow, i.e. 


n = 


( 2 ) 


In terms of the following dimensionless variables 

x = t/(d 2 /v), (x j z) = (X j Z ) /d , ip = V/v , £ = n/(v/d 2 ) , 

where d is the characteristic length, (1) can be written as 


+ V X " Vz = V « • 


(3) 


S - . 


Taylor (1923) obtained an exact solution of (3) given by 


iji = A cos ttx cos ttz exp(-2TT x) . 


(4) 


where A is a constant representing the strength of the flow. Without loss 
of generality, A is taken to be positive. This stream function represents 
a time dependent two-dimensional array of counter-rotating vorticies 
enclosed in squares of a dimensional side-length d. Note that 
x - (K+D/2 and z = (K+l)/2 are stream lines, where K is an integer . 
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Note also that the x- and z-components of the velocity are respectively 
given by 


u = \|> 2 = -Attcosttx sintfz exp(-2TT t) 


and 


v ® -\p^ = ATrsinux costtz exp(-27T x) 


Hence the above mentioned two sets of stream lines intersect at saddle 
points, and the centers of vorticies are located at the stagnation points 
(x,z) = (M,N) where M and N are any integers. 

The nonlinear stability analysis of this flow does not seem to have 
yet appeared. Most of the exact solutions of the Navier-Stokes equations 
represent parallel flows. The linear stability of parallel flows is giverned 
by the well-known Orr - Somerfeld type ordinary differential equation. 
(Joseph, 1976; Drazin and Reid, 1983) The linear stability of the 
non-parallel flow described by (4) is giverned by the partial differential 
equation 




where £ = V ip , 
and vorticity. 


- £ C 1 + ^ ' C - ' ? = V 2 ? ' 

X Z Z X x z 


(5) 


and ip ' and are respectively the perturbation stream function 
(5) is obtained by substituting 


ip = ip + w 1 and £ £ + £ 1 

into (3) and retaining only the linear terms. Based on this equation, it 
was shown recently by Lin and Tobak (1986) that Taylor's vortex array is 
absolutely and monotonically stable with respect to two-dimensional infinitesimal 
disturbances which vanish at infinity for all frequencies. 
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It turns out that Taylor's vortex array is also globally stable with respect 


to three- as well as two-dimensional disturbances of any finite amplitude. 


This will be proved presently. The Navier-Stokes equations can be written as 


3 t (U+u) + (U+u)-V(U+u) 


1 2 
~V( p +p)4vV (U+u) 


V • (U+u) = 0 , 


where U is the basic flow field the stability of which is of interest, u 


is the velocity perturbation, p is the density, p is the pressure field in 


the basic flow and p is the pressure perturbation. The basic flow field 


satisfies 


3 t U + U-VU = 


-Vp + vV^U 
P“ ~ 


v-u = 0 


Substrating (8) from (6), forming the inner product of the resulting equation 


with u and using the fact that the velocity field is divergence free, we have 


[d t + (U+u).V]e(t) = -u(u-7)U-V. (up) + { V * (u . V )u ] -7u : Vu} 


where e(t) = i s the kinetic energy of the disturbance per unit mass 


of fluid. Integrating (10) over the entire flow domain, using the divergence 


theorem and the boundary condition 


l\ s = °» 


where the subscript s denotes that u is to be evaluated at the boundary s of 


the flow domain, we have 


= -v<Vu:Vu> - <u*(u«V)U> , 


( 12 ) 
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where d/dt = 3 t + (U+u) • V . For two dimensional basic flows, such as the 
Taylor vortex array, we have 

<u-(u-V)U> = <u?3 x U + v 2 3 z V + uv(3 x V + 3 Z U)> , (13) 

where U and V are respectively the X- and Z-components of the basic flow 
field U, and u and v are the corresponding components of the disturbances. 

It should be pointed out that u was not assumed to be two-dimensional in 
(13). It can be seen from the Taylor solution (4) that the shear rate of 
deformation in the basic flow vanishes, i.e., 

9 V + 3 U = 0 , (14) 

A L* 

and 

3 X U = =v (iT/d^)Acos(7TX/d) . cos (irZ/d)exp( -2 tt^v t/d^) . 

2 2 2 

Hence the minimum of 3^U and 3^V is -(v^/d )Aexp(-2?r vt/d ). It follows 
from (13) and (14) that 

<u • (u- V )U> -< (vtt / d )(u^+v^)Aexp(-2Tr^vt/d^)> _> -(vir/d^) . 

<(u^+v^+w^)>A exp(-2Tt^v t/d^) . (15) 

On the other hand, it can be shown by use of Schwarz's inequality that (see 
e.g. Joseph 1976) 

<Vu:Vu> > <(u-u)> , (16) 

I 1 

where 1 = Kd, K being any arbitrarily large integer, is the distance 
between any arbitrarily far apart parallel boundaries of the Taylor 
vortex array where condition (11) is invoked. It follows from (12), (15) 
and (16) that 
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j- <E> < -2<E> + \ A exp(-2TT 2 vt/d 2 )] . (17) 

JT d Z 

Integration of this inequality yields 

<E(t)> < <E(0)> exp[ I( t ) ] , (18) 

I(t) = - [1 - exp(-2ii 2 vt/d^) ] , 

K 2 d 2 * 

where <E(0)> is the energy integral of any kinematically admissible initial 
disturbances. Since the inequalities (15) and (16) used to arrive at (18) 
hold for both two- and three-dimensional disturbances, the decay law 
described by (18) holds for both types of disturbances. 

In summary, we state the following results. Taylor's two-dimensional 
vortex array, which is an exact solution of the Navier-Stokes equations, is 
asymptotically stable in the mean in a sense that 

Lim <E( t ) >/<E(0 )> + 0 (19) 

£-><» 

with respect to two- as well as three-dimensional disturbances satisfying 
the kinematic conditions 

V-u = 0 , u| s = 0 , v > 0 and <E(0)> < <=° . (20) 

Moreover, the right side of (17) is negative for all t. Hence the 

disturbance energy integral decays monotonically . The stability is global 
in a sense that the initial disturbance is allowed to exist everywhere in 
the infinite flow domain without limitation except that given by (20). This 

can be seen easily from the fact that (19) remains true even when K-*». It 

can be shown easily from (4) and (18) that as K-*», the initial decay rate of 
the disturbance is greater than that of the basic flow if A>tt . If A<tt , 
the converse is true. It follows from (18) that the decay rate decreases 
as K and d increase, but never vanishes for finite d. 
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Finally, we note that this work, as far as the authors are aware, 
is the first nonlinear stability analysis of a non-parallel flow which 
is an exact solution of the Navier-Stokes equation. For this reason 
this work will probably serve as a valuable testing ground for any method 
of stability analysis of non-parallel flows for which an exact solution 
is not available. Unlike the parallel flow theory, the general method 
of stability analysis of non-parallel flow is yet to be firmly established 
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Introduction 


Most of the known exact analytical solutions of the Navier-Stokes 
equations have been obtained for parallel laminar flows for which 
the Navier-S tokes equations can be linearized. G.I. Taylor*' 
observed that the nonlinear convective terms in the two-dimensional 
Navier-Stokes equations vanish when the vorticity is a function 
of the Stokes stream, function alone. For the special case of 

vorticity being proportional to the stream function, he obtained 

\ 

an exact solution which represented a double infinite array of vorticies 

2 

decaying exponentially with time. Kovasznay extended G.I. Taylor's 
method, and considered the flow in which the local vorticity is 
proportional to the stream function perturbed by a uniform stream. 

He was able to linearize the Navier-Stokes equation and obtained 
an exact solution which may be used to represent the flow down stream 
of a two-dimensional grid. 

In this note we make a minor extension of Kovasznay' s solution, 
and present an exact solution of the Navier-Stokes equations, which 
may represent the reversed flow above a flat plate with suction. 

The present exact solution grew out of our search for a simple exact 
solution which may serve as an exact basic flow for a model stability 
analysis of non-parallel flows involving a flow reversal. Unlike 
the stability analysis of parallel flows, the stability analysis 
of non-parallel flows is not yet firmly established. Among other 
difficulties, the stability analysis of nonparallel flows suffers 
from some serious uncertainty associated with the approximate nature 
of non-parallel basic flows. 


I 



Exact Solution 


Consider two-dimensional flows of an incompressible Newtonian 
fluid. The X and Z Cartesian components of the velocity field are 
related respectively to the Stokes stream function V by 

U " T Z ’ V ’ " *X ‘ 

where subscripts X and Z denote partial differentiations. In terms 
of the governing equation is 

O -vV 2 )7 2 'P = 't' v 7 2 'f_ - '? 7 V 2 4 ' v , (1) 

2 

where t is time, v the kinematic viscosity, 7. Laplacian 

operator, the subscripts X and Z denote partial differentiations 
with respect to the space variables. By use of the following dimension- 
less variables 

^ = y/v , ( x , z ) * (X,Z)/(v/U) , 

7 = 7(v/U) , T = t/(v/U 2 ) , 

where U is the characteristic velocity, Eq. (1) can be rewritten as 

(3 - V 2 )V 2 l> * i|) v V 2 ii) - 1 1> V 2 i i> . (2) 

T X 2 Z X 

This nonlinear equation can be linearized for flows with the following 
particular vorticity distribution 

7^ = k(0-Rz), (3) 

where k is a constant to be determined, and R is a flow parameter to 
be explained. Substitution of Eq . (3) into Eq - (2) gives the following 

linear ^quac ion 




i 


( 4 ) 


(3 - V 2 )t|> = “fty • 

T X 

A special steady solution, ^ , of Eq . (4) is given by 
^ = Rz + B exp(mx-nz) , (5) 

where B and n are constants, and 

m = %(R + /R 2 -4n 2 ) . (6) 

The coefficient of vorticity distribution, k, appearing in Eq. (3), 
can now be determined by substituting Eqs. (5) and (6) into Eq . 

(3). It is found that 

k = mR. 

When n=2*n’i, and R is identified as the Reynolds number, the 
solution given by Eqs. (5) and (6) reduces to the Kovasznay flow. 

Kovasznay suggested that the solution corresponding to the negative 

.... ... . — 

2 2 

root of m in Eq. (6), i.e. ra = 0.5 (R-/R +16tt ), can be used to 

describe the wake flow in the region x>0 behind a grid at x=0. 

The solution corresponding to the positive root of m, i.e. 

2 2 

ra 3 0.5 (R+/r +16tt ), can be used to describe the two-dimensional 
uniform flow, in x<0, perturbed by a grid of blunt bodies with suction 
or blowing whose magnitude is proportional to B near the front stagnation 
points . 

When n is real, the new exact solution given by Eqs. (5) and 

(6) represents flows over a plate with suction or blowing of fluid 

2 2 

at the plate. For example, when n>0 and R -4n >0 the exact solution 
describes the flow over a porous plate z=0, x<0 with velocity components 


* 


3 



given by 


ORIGINAL .PAGE' IS 

u = ip - R ~ nB exp(rax-nz) , OF POOR QUALITY (7) 

z 

v = —((i = -raB exp(mx-nz) . (8) 

It is seen from Eqs. (7) and (8) Chat R can be put to 1 without 

loss of generality, since it amounts to rescaling the velocity by 

a factor of R. It is seen from Eq. (8) that B>0 corresponds to 

suction and B<0 corresponds to blowing at the bottom plate z*0. 

Fig. 1 gives the streamlines of the flow for R=l, n=0.05 and B=80. 

Note the flow reversal due to suction near the trailing edge. It 

should be pointed out that the flow reversal may not take place 

when B is sufficiently small. 

2 2 2 

When n>0 but R -4n = -q <0, the real part of Eq. (5) gives 

* Rz + B exp(Rx/2-nz) cos qx. 

Hence 

u = R-nB exp(Rx/2-nz) cos qx , 

v = -B[(R/2)cos qx - q sin qx] exp(Rx/2-nz) 

This solution represents the wavy flow over a plate z=0, x<0 with a 
periodic variation of suction and blowing along the place. The possible 
flow represented by Eqs. (5) and (6), with n neal, are delineated in Fig. 2 
Note that R is put to 1 without loss of generality. 

Discussion 

A new exact solution of the Navier-Stokes equation is given. 

The solution can be used to describe the flow with or without flow 

reversal above a plate with suction, and/or blowing. Aside from 
its own intrinsic interest, it provides us with a simple but exact 



non-parallel basic flow which could be used for a model analysis 


of stability of non-parallel flow. The stability analysis of non- 
parallel flows usually suffers from the uncertainty associated with 
modeling the non-parallel flows with approximate solutions. 
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In l roduc r. ion 


\ . 

The laminar boundary- layer oscillation on a rotating disk was 
investigated by Smith (1947). The spiral vortices on a rotating 
disk were discovered by Gregory, Stuart & Walker (1955). They 
theorized that the sprial vorticies were the products of flow 
instability. Their study has been expanded and improved by many authors 
Gregory & Walker, I960; Brown, 1961; Tqbak, 1973; Kobayashi, 1980). The 
agreement between theories and experiments remain inconclusive. Similar 
spiral vortex motions are observed on other surfaces where the centrifugal 
force plays an important role; spiral vorticies on a swept wing are an 
example (see Chapman & Tobak, 1984). The boundary-layer and side force 
characteristics of a spinning axisymmetric body were investigated 
experimentally by Kegelraan, Nelson, & Muller (1980) and Muller et al. 
(1981). They observed, in addition to the Tollmien Schlichting waves, 
regular vortices roping tightly around the rotating cone-cylinder. 

Similar vortex-motion has been observed by Kobayashi (1980), Kobayashi 
and Izumi (1983), Kobayashi et al. (1983) on a rotating cone. They also 
analyzed the flow with the aid of linear stability theories. Definitive 
comparisons between theories and experiments cannot yet be” made. 

The lack of definitive comparisons between experiments and theories 
mentioned in the previous paragraph is not surprising. This is because 
in all of the theoretical studies, flow instability was assumed to be a 
local phenomena. The basic flows were "frozen" locally in the analysis 
which disregards the possible up-stream and down-stream influence. The 
concern over this type of analysis for non-parallel flows has already 
been raised by Morkovin (1979). Note that all of the basic flows 
mentioned in the previous paragraph are non-parallel. The stability 


theories Co; parallel t Lows am far more well developed than the theories 
for non-parallel flows (Joseph, 1976). The parallel basic flows possess 
translational invariance in the direction of flow. No such invariance 
exists for non-parallel basic flows. While it is reasonable to assume 
that the physical characteristics of the disturbance retains the same 
invariance in the case of parallel flows, no such assumption can be made 
for the case of non-parallel flows. The mathematical consequence is that 
the governing equation of stability is an ordinary differential equation 
for the former case but it is a partial differential equation for the 
latter one. The analysis of the latter case is made even more difficult 
by the fact that while the parallel basic flow is independent of the 
parameter, the critical value of which is to be determined by stability 
analysis, the non-parallel basic flow depends crucially on the stability 
parameter itself. Stability analysis for non-parallel flows is still in 
an exploration stage (Goldstein, 1981; Goldstein et al., 1983; Hall, 1983). 

The purpose of this work is to investigate the precise mechanism of 
the flow separation from a rotating body of revolution. The spiral vortices 
around a rotating body are viewed as the consequence of the instability 
of the steady non-parallel laminar basic flow about a body of revolution. 



Basic Flows 
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Two different basic flows are considered. The first is the fLow 
around an axisymmetric body rotating at a constant angular velocity 
ft in an otherwise quiescent fluid of infinite extent as shown in figure 
1. The second is the flow around the same body placed in a coaxial 
cylinder of infinite length. The axial stream far from the nose-tip 
may be imposed externally for this latter flow. For sufficiently small 
Reynolds numbers based on Qa, a being the maximum radius of the body, 
the flow is laminar and axisymmetric. For these basic flows the Navier- 
Stokes equations are reduced to 


DU _ vf 
DC ’ r 


9 P 
r 


l 

R 


(v 2 u - V 


DV UV 
Dt r 



DW 

Dt 


a p + 

z 



J. 

R 


v 2 w 


* 


3 (rU) + 3 (rW) = 0 , 
r z 


( 1 ) 


( 2 ) 

(3) 

(4) 


where t is the time nondimensionalized by Q \ (U,V,W) are the radial, 

tangential and axial velocity components nondimensionalized by fta 

in Che cylindrical coordinates (r,0 t z); all lengths are normalized with 

2 

a, R is the Reynolds number Da /v , v being the fluid kinematic viscosity, 
and 


D_ 

Dt 


3 + U3 + - 3 + W3 

c r r 9 2 


(5) 


7 


2 


3 

rr 




96 


+ 3 

zz 


(6) 


The pressure terms in (1) and (3) can be eliminated by cros s -d i f f erent iat ion 
which yields, upon utilizing (4), the vorticity equation 



( 7 ) 


HL 

Dt 



2 V 1 

— 3 v + i (v* 

r z R 



r 


where £ is Che vorticity given by 


£ = 3 U - 9 W . 

z r 


The velocity components can be obtained from (2), (4) and (7). Then* 
if desired, the pressure field is obtainable from (1) and (3). 

A Galerkin method similar to that used by Munson & Joseph (1971), 
Munson 6c Menguturk (1975), Yakushin (1969, 1972, 1973) and Khlebutin 
(1968) will be used to construct the solution. The velocity components 
will be represented by a complete set of orthogonal functions which 
satisfy the no-slip boundary conditions on the rigid surface* For 
the convenience of applying the boundary condition, we introduce the new 
variable y 3 r-h(z,x). Thus r and 3 r in (1) to (6) must be replaced 
by y+h and 3y respectively, and the differentiation with respect to 
z must be transformed according to 


(3 ) - D - h'd 

z r ,x 

2 

O ) » (D-h'd)(D-h'd) = D 2 -2h'dD + h' d 2 - h"d, ecc. 

zz r ,x 

where primes denote differentiation with respect to z, and 

D - ( 3 ) and d = (3 ) 

z y.T y z,r 

2.1. Unbounded domain 

Consider the steady axisymraetric flow around a body of revolution 
in an unbounded fluid as shown in figure 1. The veLocity components are 
represented by a complete orthogonal set 


V 


■' ! Vi.i, = ' GftlGSNAL PAGE 2S 
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J 

V - F = . VnUJ Yn(y) , F = h‘/( v+h> (8) 


W 


= Wn( z ) Yn ( y ) , 


where the double indecies imply summation over n, from n-1 to N, and 
f s -s 2, s f v 

(ye ) L (y) 

V y) • uTr(7^fT7^ ■ s> -' • 


(9) 


(y) in (8) is the generalized Laguerre polynomial ( Aizenshtadt , 1966) 

possessing the following orh cogona 1 i ty relation 

00 - e 0 , n f m 

/ y S e 3 L S L dy = { 
o n l P(s+n+l ) , n = m ; 

and T ( s+n+ 1 ) denotes the Gama function. 

f ( s+n- 1 ) = ( s+n ) ! . 


Note that the’ velocity components given in (7) satisfies the boundary 
condition at Che body surface, 


U = W = 0 and V = h at y ■ 58 0 , 

and the condition that the velocity vanishes at infinity. 

Substituting (8) into (2), (4) and t ,7 ) , we have 

2 

{ ( ( D-h ' d ) “+d +d/ r-r* “ ] -R( U_ Yd+W Yn(D-h'd)+U Y„/r]}(V Y +F) = 0 , 

n n n n n n n 

dU Y +U Y n / r+( D-h ' d )W V = 0 , 
n n n n n n 

2 

{ ( ( D-h ' d ) 2 +d +d/r-r' 2 )-R[U Y n d-fW Y (D-h'd)-U Y /r J } ( (D-h ’ d )U_ Y n -dW y 1 

n n nn nn nnnn 

+ ( 2R( V y n +F )/r](D-h'd)(V Y +F ) = 0 . 
n 11 n n 


* - 
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Multiplying ach of the previous three equations with Y , Integrating 

01 

over 0 < y < °°, and using Che orthogonality property of Y we have 
— ' " m 

V "6 =2h'V ' / Y Y ' - V / Y d 2 Y 

n mn n m n nmn 



+ R ft Y [U dF + W (DF - h'dF) + FU /r] 
n .m n n n 

- / (d 2 +D 2 -2h'dD)FY 

m 


V 


6 = h ' W n /Yv 1 

mn . n ® n 


U[/Y Y '+/(Y m Y /r)] 


nmn 


U M, 6 
n mn 


U "*3h'/Y Y ' - U '[2h' 2 /Y Y "-2h"/Y Y '+/Y d 2 Yn] 
” m n m n m 


n 


m n n 


U [-h"'/Y Y dy-h'/Y d 2 Y +2h'h"/Y.Y "] + W " f Y Y ' 
n mn mn J mn n'mn 


W_ ' *2h'/Y Y n " + W/Y d^Y 


2 V , 


m 


n m n 


( 10 ) 

( 11 ) 


+ R{U "-W Y Y Y_ 
n e m e n 


U * [U / Y Y_Y n ' 
n e me n 


- 2 W h'/Y Y e Y 
e m • n 


- U /Y Y Y /r 
e m e n 


- h 1 W /Y Y Y '] 
e ta e n 
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+ U[.-Uh'jY Y e Y "-h'U fr" 1 Y Y e Y '+h' 2 W Jy Y e Y " 
n l e J m e n . m c n e J m c n 


W h"/Y y Y n ' ] 

e nj 1 e 11 


+ M n ’[-W e /Y B V e Y„') 


W [ U /Y y Y "*U /r -1 Y "Y Y 1 -h'U./y Y.Y "1 
n e m A e n e men e A m e n i 


- 2V * [V /Y Y +/r FY Y ] 
n e me m n 


+ 2V [h’V /Y YeY , +h’/r~ 1 FY Y_' 
n e m n m n 


- f Y Y DF+/Y Y h'dF] 
m n m n 


- /r' 1 Y F(DF-h'dF) = 0 , 

m 


( 12 ) 


where primes denote differentiation with respect to the arguments of the 
functions, the integration is understood to be between y=0 and each of e 

m and n range from 0 to N, and 


d 2 


+ - + h' 2 d 2 - 


h"d - r - 


(10), (11) and (12) are a system of 3N equations in 3N unknowns Un, Yn 

and Wn (n a l to N). The independent variable is z. This system is third 

order in U , second order in V and first order in W . NoTe the second 
n n n 

order term of W in (12) can be reduced to the first order terms by 
n J 

use of (11). This system can be reduced to a system of the first order 
differential equations given by 



iV 15 'J 

m la 

V l / 

t 



! I ) 

(u ' lJ r 

m 

- u (2) , 

m 

• 


(u (2) r 

m 

- s (u , u < l >, u (2) , V V (l) , W , z ) , 

Inn n n n n 


CL3) 


(L) 


V ' = V 
m ra 


[V (1) 1' = s,(u , V , v (l) , W , z) , 
m Z n n n n 


[W J 1 = S,(U , W f z) . 
® j n n 


where S^, S 2 and are given respectively by the right sides of (12), 

(10) and (11). The solution of this system with appropriate upstream 
conditions has been attempted by the use of the Runge-Kutta method, for 
various values of R. Seveor numerical instability has been encountered, 
possibly because of the extreme spatial resolution required. Note that 
(13) is a dynamical system in the variable z. The fixed points of the system 
in the phase space are given by the solution of (13) with all elements on 
the left side put to zero. Note that this system is non-autonomous , since 
z appears explicitely through the body surface r - h(z). Thus the character 
of the fixed points is expected to depend on z. To overcome the numerical 
instability encountered with the Runge-Kutta method, we expand the amplitude 
functions in (8) in terms of another complete orthonormal set Pk(z) and 
apply the Galerkin method in the z-direction. Hence 


(U(z), V (z), W (z)] = [U , , V W ] P (z) 
n n n nk nk nk k 


(14) 


where the coefficients of P^(z) are constants, and the repeated subscript 
k denotes summation over k. Substituting (14) into (10) to (12), multiplying 
the result ing equations with P ( z ) , and integrating from z 3 -« to z=+®, we have 


V 



i5 

rnn 

6 

mn 

6 

mn 


I V . P, "P, = j Si (m , n . k ) !* 
j uk k ?. J 1 l 

j w i.*\ = I S, (m, n , k ) P 

; mk k l J 2 I 

j u P *"P = J S 3 (m,n,k)P. 

1 nk ^ J ’ * i 
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( 15 ) 


where S^, and are the respectively the right side of (10), (11) and (12) 

with the velocity amplitude functions properly expanded in terms of P •. 

can be any set of orthonormal functions which vanish at infinity and render 

the right side of (14) integrable over -«<z<°°. (15) is a system of nonlinear 

algebra equations in the amplitude functions U, , V and W, . The 

kn kn kn 

solution of this system can be achieved in principle by the Newton-Raphson 
method. For the convergence of this method, as usual, a sufficiently close 
initial guess of the solution is required. For steady nonlinear flows, 
such a close initial guess is difficult. To avoid this difficulty, we may 
treat the amplitude functions as time dependent and solve the initial value 
problem. The solution then may be started from any kinematically compatible 
initial condition without the need of a close initial guess. For the 
initial value problem, we need only to add to the left sides of the first 
and the last of (15) respectively the following terms, 


-<5 { u, 3 » 

mn ki t kn 


and 


I I V W 


nk * 


The resulting equations form a system of ordinary differential equations in 

[U (t ) , V (t ) > W (t ) 1 . The solution of this system can be achieved by 
nk nk nk 7 

use of the finite difference method successfully employed by Aidun (1985). 


2.2 


3ounded domain 


For the tiow around a bodv of revolution, bounded intcrnallv by a 

V. ' 

cylinder of radius b as shown in figure 2, we represent each velocity 
component by an orthogonal set, each member of which satisfies the no- 
slip condition at solid walls. Thus 


U = U ( z)Y (y,h) 
n n 

V-F “ V (z)Y (y,h> (16) 

n n 

W-G = W (z)Y (y,h) 
n n 

where 

Y n = sin nn 

F = A(y+h) + B/(y+h) , A = -h 2 /(b 2 -h 2 ), B = b 2 h 2 / (b 2 -h 2 ) ; 

G = ( 1-H) ( l-y/b) EnTrVM-®) 

1 . , z>0 

H ■ ?'' 

0 , z<0 


Note that G is so chosen that the axial velocity is allowed to exist at 
z=-H» , and dW/dy = 0 in z £ 0. F is so chosen that the circumferential 
velocity at z-« approaches chat between two concentric cylinders with 
the inner cylinder rotating and the outer one fixed. The forms of Y^ and 
F allow the no-slip condition to be satisfied at the solid surfaces. 

Substituting (16) into (2), (4) and (7) and integrating between y=0 
and y-b-h with the weighting function Y 

m 

V "6 ( b-h ) / 2 = V 1 ( 2h ' fdY Y -2(DYY ) 

n mn n ; n m 1 n m 




( 17 ) 


- h ' V? /Y ,dY Y -h'/GdY Y-+/GDY Y m 


+ U /r _1 Y Y Y ) + [W /Y Y (D-h'd)F 
e e n m n n m 


+ /CY (D-h'd)F + U /r' l FY n Y +U / Y Y dF} 


m 


ra n n 


W '6 (b-h)/2 = W [H'/dY n Y +/DY_Y m J-U n ((dY n Y +|r' l fY ) 

nmn n J n m n l ni n; nm J nra 


- /(D-h')GY 


m 


U "'5 (-bh)/2 * 3U " [h ' JdY Y +/DY Y ]+3U ' [ -/ (d 2 )Y Y 

n rnn n n m n ra n n m 


3. 


- K Y n Ym + 2 h ’^DY n Y m ]+U n (-/d 1 Y^ 


/D 3 Y Y +3h'/dD 2 Y Y -3/(d 2 )DY n Y 


n *ra 


n m 


W '7dY Y +W l-2h’/d 2 Y Y m +2/dDY Y m ] 
n n m n n 


+ W [/d- 2 dY Y m +/D 2 dY Y--2h'/d 2 DY Y 1 
n 2 n m n m nm 


/(D-2h'dD+d 2 )dGY 


m 


+ R{-U "[W /Y Y Y +/GY Yl 
n e e n m n ® 


+ O' (2h'W e /Y e dY Y +2h'/GdY Y-+U /Y dY Y 
n e c n ® a“e e nm 


-1 


V r Y e Y nY m -2W e /Y e DY n Y m -2/GDY n Y m ] 


+ U [U / Y dDY Y m -U /r' L Y DY Y -h'U /Y d *Y Y 
n ee n“ e e nm e e n® 


-1 


-1, 


+ h'U /r Y dY Y +/Y n d GY -/r Y dGY 


e e n m 


m 


n m 


+ 2h 'W /Y DdY Ym+ 2h 1 /CDdY Y 
e e a 111 n m 


- W /Y D 2 Y Y m -/GD 2 Y„Y 
e e a ra am 


+ W ' (W /Y dY Y +/CdY Y ] 
n e e n m n m 


(18) 


where 


-l 


+ W (U fr Y dY n Y m -/Y e d i Y +W / Y OdY Y 
n e enm e n m e e R m 


+ /GDdY n Y m -W e /Y e h ' d 2 Y n Y m -/Gh ' ( D ' h ' d )dGY m 1 


-I 


- 2V ' [V /r *Y Y Y +/r~ l FY Y J-2V [V / r ' 1 Y DY Y +/r' L FDY n Y 
11 enm n m ne e rt m n 


-i. 


- h ' v e ^ r " 1 Y e <iY n Y nI -h , /r"^FdY n Y aI 

+ /r’ l Y (D-h'd)FY *2/r’ L F(d-h'd)FY } 
am m 


(19) 


5 . 3 


L 3 * -(h’+h' 3 )d 3 +( 3 h , h"-r" l h' )d 2 -h’"d 


Again the differential system can be casted into a system of 6 N nonlinear 
ordinary differential equations of the first order of the same form as 
(13) where and S 3 are now respectively given by the right aides 

of (19), (17) and (18). 

The extremely high spatial resolution required for the Runge-Kutta 
method to converge was also encountered in the bounded domain. The 
alternate methods of minimizing the residual over the entire flow domain 
and treating the flow establishment as an initial value problem described 
for the unbounded domain can be equally applied to the bounded case. 


3. Stability Analysis 

The solution to the differencial system given in the previous subsection 

depends on the flow parameter R. The assumed ax i symmetric flow may not 

exist due to instability at Reynolds numbers exceeding a critical value 

R . The R at the onset of instability may be, in principle, determined 
u c 

from extensive repetitive solutions of the differential system for different 
increasingly Large R values until the ax i s ymrae t r ic solution becomes 




I 


impossible. New bifurcated non-ax i a ymrnet r ic solutions must be found. 
However, such computation is costly and not physically very revealing. 
For this reason, we will investigate the vortex structure in the 
separated flow as the consequence of the instability of the basic flow 
given by ( 8 ) and (>16). 

According to the known observations (Muller et al., 1981, Kobayashi 
et al., 1983), the perturbed flow is not axisynnnetric. Hence the 
governing equations of the non-axisymme trie flow perturbations are 


Du 2Vv i 7 I, 

TT* + U + v3 U - — — * -3 P + T [v u - 

Dt y z y+h y ^ R , v 2 


( y+h ) ( y+h ) 


2 0 


3* v l * 


Dv ^ ^ \ f , , Uv+Vu 

— + u3 V + w3 V + — 

Dt y z y+h 


y+h 


3 qP + 


1 f7 2 

- |V v 


v 2 

+ =■ 3 0 u] , 


2 ' 2 9 

( y+h ) ( y+h ) 


^ + u 3 y W + w3 z W = 


3 p + 7^w 

z r R 


C 20) 


3r(ru) + 3 (rw) + 3 v 13 0 

■z : 9 


where (u,v,w) are perturbations of (U,V,W). The continuous septrura of the 
above system is assumed to be damped, and the discrete spectrum which 
portrays the spiral vorticies is periodic in 9. Hence, we represent 
the solution by the complete set 

u = u ( r , z ) Y exp [ 180 ] , 

n n 


v = v ( T , z)Y exp( 169 1 , 

n 


(21) 


W * w n (T^ 2 )Y n exp[ 180 ] , 


where 8 is the real wave number in the tangential direction. Note that 
( 21 ) satisfies the no-slip condition at the solid surface and/or the 
vanishing of disturbances at y*«* . Substitution of ( *J i ) into (2U) and 



e 1 ituina t ion of the pressure terms by cross-dif f eren t iat ion then yield a 
system of 3N partial differential equations with (U^, v n » W^) as 
dependent variables and x and z as independent variables. Solution 
of this system with an appropriate initial condition and an upstream 
condition but without specifying the downstream condition yields the 
information sought about the instability. 

An alternate approach is to solve (21) as an eigenvalue problem. 

To do this, we expand the amplitude functions in a complete orthonormal 
set P^(z), 

u 3 u P (z)exp(ojT ) , 

n nk k 

v 3 v P (z)exp(cjt ) , 

n nk k r 

w = w P (z)exp(orc) , (22) 

n nk k 


where u> is the complex eigen- frequency . Substituting (21) and (22) into (20), 
multiplying the resulting equations by P^ and integrating from z a -<» to +°®, 
we obtain a system of homogeneous linear equations in [u^> v n ic J W nk^ * 
existence condition of the solution then determined the eigen- frequency 
which determines the stability criteria. The details of this approach 
are demonstrated in the accompanying report "Stability of the Kovasznay Flow." 
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I. INTRODUCTION 


THE unsteady flow separation about a radially deforming cylinder 

starting an arbitrary translation from rest is analyzed. This 

problem is relevant to several fundamental issues of considerable 

importance. These issues include the control of flow separation 

by moving an aerodynamic surface^ - the possible drag reduction 

2 

by surface deformation and the issue on the appropriate 

3 

definition of an unsteady flow separation point . 

However, this study was originally motivated by 
the need to understand the origin of flow separation in 
a steady three-dimensional flow about a body of revolution translating 
with an angle of attack in an otherwise quiescent Newtonian fluid. As 
one moves downstream along the axis of revolution and views the flow 
pattern on successive planes perpendicular to the axis, one sees the 
evolution of a sectionally two-dimensional vortical flow structure 
about a circle whose radius varies with time. Thus there is a 
certain analogy between the steady three-dimensional flow about a 
slender body of revolution and the unsteady two dimensional flow about 
a circular cylinder whose radius changes with time in a cross-flow. 

Uf, 

This topological analogy was demonstrated by Tobak and Peak , 
and Chapman and Tobak^. 

Creeping flow around a deforming sphere was studied by 
2 

Lin and Gautesen on the basis of the unsteady 

Stokes equation. The linear unsteady Stokes equation is known to 

be also a valid first order approximate equation even at finite Reynolds 

6-8 

numbers for small amplitude high frequency fluctuating flows 

9 

as well as for initial transient flows . The nonlinear convective acceleratic 


i 



terms were later treated as source terms in perturbation solutions to 

the Navier-Stokes equation by Tseng and Lin^’^ in their 

study of transient heat transfer from a heated wire and for constructing 

a theory of a heat sensing velocimeter. While these perturbation 

solutions are useful for describing the initial transience of flows in 

which the convective acceleration is weaker than the local acceleration, 

they are obviously inadequate for the description of separated flows 

in which the convective and local accelerations may become equally 

important at a later stage of the flow development. 

An iterative scheme of solving the Navie,r-Stokes equations for 

unsteady two-dimensional, flow about a deforming and translating 

cylinder is given in the next section. The method is novel to the 

12 

best knowledge of the writers . In the first 

approximation, the effect of vorticity convection is neglected. This 

effect is taken into account by successive iterations. Thus, in the 

k-th iteration, the convected vorticity appears as the source term 

in an equation of transient diffusion of vorticity. A novel integral 

transform is used to reduce this transient vorticity equation to a 

k-dimensional heat equation. The bounded solution of this equation with 

the appropriate boundary conditions is given in section 3. A general 

method of superposition for problems involving a moving boundary is 

given in the Appendix. An equation which describes the migration of 

the separation point on a deforming cylinder is obtained in section 

4. The numerical results are then presented and discussed. The present 

iterative solution is expected to be valid during the initial transience 

before a large amount of vorticity is convected to infinity. Then the 

13 

well known difficulty associated with the Whithead paradox 

surfaces. The numerical results show that the flow separation establishes 
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itself long before the vorticity has time to populate itself at infinity. 


II . Formulation 

Consider the flow of an incompressible Newtonian fluid about a 
circular cylinder. The fluid is initially quiescent. The translational 
velocity of the cylinder is -iU(t) , where t is time and i^ is the unit 
vector in the negative x-direction as shown in Fig. 1. The instantaneous 
radius of the deforming cylinder is given by a^(t). 

The governing equations with respect to a reference frame attached 
to the cylinder are 

7*V = 0, 

(1) 

3 t I + Til — j vp + w 2 v + i u(t) - j.g , • 

where V is the velocity, P is the pressure, p is density, v is the 

kinematic viscosity, and j_g is the gravitational force per unit mass of 

fluid. Let the characteristic time and distance be respectively ui ^ and 

6. The corresponding characteristic velocity is . Nondimens ionalizing 

2 

the distance, time, velocity and pressure respectively by a g , a g /v, 

2 

u)6 and p(w6) , Eq. (1) can be written as 


3 v 

T— 

+ eR(v* 7) v = 

-eR7p + 7 2 v 

9 


(2a) 

7* v 

= 0 

9 




(2b) 

where 







v = 

V/o5 , 

T « 

t/[a o 2 /v], 

e 

* «/>o- 


P = 

(P + 

gy - u(t 

:)x) /p (co6 ) 2 , 

R 

■ 2 

- 0)a o /v , 


and a i.- 
0 

; the 

avern.ue 

cy ] inder rad 

ins 

i . e . , a = (a . + a ) / 1 . 

o mm max 



Equation (2b) is the necessarv and sufficient condition for the 



exiscence of the stream function ip such that the radial and tangential 
velocity components of two-dimensional velocity can be written respectively 
as 


v = 3 'b , 

r r 0^ 


v. = 3 <jj 
0 r 


In terms of the stream functions, Eqs. (2a) and (2b) can be written as 


O x - V 2 )v 2 ^ 


R 3(»,7 j) = _ R 
r 3 (r,9) £ r 


(3) 


where 


2 11 
7 * 3 + - 3 + d aa , 

rr r r 2 68 

r 

J = 3 r 03g7 2 0 - 3 0 O3 r 7 2 O • 


The initial condition is 


0(r ,8 ,t) = 0 , t < 0 


The boundary conditions are 


- ( 7 V’r-aCr) 'VWt) * 0 1 


and 


" ( r 3 0 4 ' ) r^o * U < T > cos 9 ’ 


-u(t) sin 9 , 


(4) 


(5) 


where u(t) => U(t)/tU(t)l 

max. 

The iterative solution of Eqs. (3), (4.) and (5) will be obtained 
for the case of symmetric flow. The generalization of the method 
to the case of asymmetric flows will be mentioned later. The 
m-th iterated solution will be written as 


5 


( 6 ) 


m 


» £ (eR) k 1 « <m) sin kO - aafl 

k-1 k 


Substituting Eq. (6) into the Jacobian in Eq . (3), we have 


(m-i) = V 2 !j/ (m 1} )/3(r,9) 


= Z I (sR)^ q+i! '" 2) {' sin q9 cos £9 
q-1 1-1 1 qZ 

m-1 


-B ( ® cos q9 sin £e}+ Z (eR) ? ~ 1 C (ra-1) sin £6, (7) 

qi J £=1 1 


where 


A ( r L) = ( 3 r ^ q (m ~ 1) )(^ 2 4-^ m 1} ) , C^ m 1} = aa3 r D^ m_1) 


ql 


(m-1) _ . (m-1) 2 (m-1). , . 

B q£ " (q *q )(a r D £^ ®- 2 

n 2 - a , -1 , .2 "2 .(0) _ n 

D„«3 + r 3 - £ r . J =0. 

£ rr r * 

By use of the relation sin q9 cos £9 = [sin (q+£)9 + sin (q-£)9]/2, we can 
rewrite Eq. (7) as 

jOn-l) = (eR) Z ~ 1 C^ m ~ 1 ‘ ) sin £9 

£-1 1 


+ H \ Z (eR) (q+t ~ 2) {(A^ _1) - B ( “" 1) )sin(q+£)6 
q-1 £=1 L qil q * 




( 8 ) 


Equation (8) can be rearranged into a form of truncated Fourier series 


b 


\ 



J (m l) = £ (eR) k-1 L) sin k« 

k-1 


+ 4 IV 3! <ER>k " 2 <a <" : ^ ■ B '” : « i)sin ke ■ hku+ "> 


l Z(m) 


, , m " 2 V I-k , _*2i+k-2,.(m-l) _ (m-1) . , 

+ ^. Z . A R A (k+4) I + B (k+4)4 

k=l 1=1 


m-2 m-1 k 2 q+k-2,.(m-l) .(m-1) * , . fl> 

+ if E I (eR) (A q(q+k) B q(q+k) } Sin ( k9) * 

k=l q=i 


(9) 


Equation (9) can be further reduced to 


J< ”’ 1) \i ca)k ' 2 { (ER)Ck "’ 1>hk ('»)' ,Js 3 1 1 [(i ' 6w)<a (” : ‘h" s ^« i> 


where 


. “’r“ k rA^- 1 ) + B (m_1) i h 

+ Z \€a>) (A / ni1 _\„ + "/o_i_i_\n7 “ 

4=1 


m-l-k 

z 


(4+k) £ ( 4+k ) l k(m-l) 

24 , . (m-1) _(m-l) 

A / A < 1 \ ** ^ rt / A I 1 ' 


t ml < £ *> + 


k9 


2 m- 2 m-1 f __ n (_.u 

+ H S S £ eR) (A (k-4)4 “ B (k-4)4 )h (k-4)(m) Sin k9 ’ < 10) 

k=m+l 4=l+k-m 


, ,1 if k<m 

\(m) ” l 0 if k>m. 


We assume that the coefficients of sin k8 with k>m are negligibly small 
in the (m-l)-th iteration, and write J^ m ^ as 


m 


j ( ~-i) . r 

k-1 


(eR) k 2 J^ m ^ sin k8 


(ID 


where J^ m ^ is given by the sum of terms in the curly bracket in Eq . (10). 
Thus, the last double summation term in Eq. (10) will be neglected. 

Substituting Eq . (6) and Eq . (11) respectively into the left and right 

sides ol: Kq . (3), and equating the coefficients of the same harmonics on 

on both sides we have 


7 



(k>l) 


( 12 ) 


2 , 2 ( m ) (m-1) 

(3 t ' Wk ' ‘ r J k 


By use of the following transform 


= u ^ T ^ a2r 1 ~ r ^ )6 i c l + ^ k 


r r 


(m) . s 2k- 1 

(S)i)s ds , 


( 13 ) 


Eq. (12) can be reduced to the 2k-diraensional heat equation 


(3 t - 7 k ) 4” ) * G k m_1> 


(14) 


where 


.(m-1) 


r r 


a(x) 


,-k j(®“ 1 )( s>T ) ds , 


7 2 = r l-2k ( r 2k *l 3 ) . 
k r r 

Note that Eq. (6) with Eq. (13) satisfies the first condition of Eq. (4). 
In order to satisfy the second condition of Eq. (4), we must have 


(a) 


1-k, 


[a( T -),t] = 2u(t )a *(r)« kl . 


(15) 


Substituting the bounded solution of Eq. (14) with Eq. (15) into the boundary 
condition Eq. (5), we have 

= u(r)cos 9 + f(r)cos k9 , 

r do 


<it) 

9r 


-u(r)sin 9 + g(x)sin k9 , 




where f(r) and g(t) are the "penalty functions” given by 


rr \ t • / -k-1 -T (m), >. 2k-l v 

f(T) = lim(r / (s,x)s ds), 

a k 

r-+°° 

g(i) = iim ( -kr / x (s,r )s ds + x, (r,T)r ) . 

a k k 

Thus the bounded soLution of Eq . (14) with Eq . (15) will not satisfy the bounds' 

condition (5) unless f(i)-K) and g ( r ) -0 . For the special case of harmonic 


8 



oscillation with m=k=l, it was shown by Lin that f(r)-H) and 

g(i)-K). In general, it is expected, on a physical ground, that there exists 

a small time x below which f(x) and g(x ) remain much smaller than u(t). In 
c 

any event this expectation can be verified -rigorously aposteriori. 

3, Solution l11 - Solution 

The solution of Eq* (14) with its boundary condition Eq. (15) can be 
achieved as follows. First, the solution of Eq. (14) which satisfies the 
quasi steady boundary condition corresponding to Eq. (15) will be obtained. 
This solution can be obtained from 


(3 T ' 7 J ) X^ m> (r,Z,T) = , (16) 

X^UUM.t] * 2uU)a 1_k (ji)<S kl H<T) , (17) 

where H(t) is the Heaviside unit step function, i.e. 

1, T > 0 , 

li(T) ■ <o, T < o . 

It is shown in the Appendix that the solution of Eq. (14) with the unsteady 


boundary condition is then given by 

x£ ra \r,T) = / 3 T x£ m ^ [r,A,T-a]d4 . 


(18) 


Equation (14) can be solved by use of the Laplace transform. The Laplace 
transform of Eq. (14) is 


[v 2 - d r 2 - (2k-l)r" 1 d r ]x 1 [ m) - g£“ _1) , (1 

2 

where v is the Laplace transform variable. By use of the new variables 
“(®) - c (m) , k-L , 

X k * /r and ? = vr , (2 

Eq . (19) can be reduced to the non-homogeneous modified Bessel equation of 
the (k-i) - th order. 
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( 21 ) 


id 2 + C^d - [l+( : k-l) 2 r~ 2 ]}S^ m) = -n (ra ~ l) r k_1 v 2 
C t k k 


The solution of this equation is given by 


^ = C b £"-< v > + B^ m) (vr)]I k _i (vr) 
+ [c< m) (v) + C^ ) (vr)lK k _ 1 (vr) 


( 22 ) 


where 


,(m) 


= -f s k K,_ 1 (vs)G^ m ^(s,v)ds, 
J 00 v 


(m) 


s k I k _ I (vs)G£ m ~ 1 ^ (s,v)ds. 


(ra) 


For the bounded solution, we chose b^ (v) = 0. Taking the inverse 
Laplace transform of Eq. (20), and applying the convolution theorem, we have 
from Eq- (20) and Eq- (22) 




(m) 


U,A) E k (r,T-A)dA 


(23) 


^03 

V, (r ,s ,t-A) G. (m “ 1) (s,A)s k dsdA} 


where q^™^ is the inverse Laplace transform of v 2 K^_^[a(£)v]c^ m ^ (v) , i.e. 


(m) 


A -1 2 


(m) 


q k = aL v K k-i [a(£)vJc k (v) > and 


V r > T) 


\ u ' x) 


£-1 !Wvr) 

^laUJvl v 


- » H t (r,T)H(T' 


-1+i 

7T 


exp [— u 


0 J, (ur)Y. - [a(A)u] -Y , (ur)J t . [a(JL)u] , 

2 , k- 1 k-1 k- 1 k-1 du 

J k-L [aU)ul + Y k-l [aU)ul U 


and is the inverse Laplace transform of 1^ ^(vr)K, ^(vs) and L^^Cvs) 

(vr) , i.e.. 
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V k = I k _ L (rs/ 2 x)exp[-(r 2 +s 2 )/ 4 T] H(t ) . 

In order to satisfy the boundary condition (17), q^ 01 ^ in Eq . (23) must satisfy 
the following integral equation 

J o (a(i) , A} H^taU) ,i-A]dA = 2u(£)S^ 

J 00 

- J o / a(jl ) V k (a(i),s,T-Al G ‘“' ;)(s ,Jl,A) sk dsdA • (24) 

» 

Note the lower integration limit in the last integral in Eq. (24) is changed 
from 0 to a(£), since vanishes for s^a. 

Recall that 


G k m L) (r,*,T) =» -J 


-k T (m-l)/ x. 

»<£> S J k U - T)<iS - 


The solution of Eq. (24r) for can be easily obtained by differentiating 

both sides of-Eq. (24) with t 

oo 

f 3 T q k ' ) [a(i),A]H k [a(Jl),T-\] [H(A)-H(A--)JdA = 2u(i)<S(t)6 kl 

•oo 

00 ® 

- 3 / [H(A)-H(X-t)]7 V (a(Z) ,s,T-A]G^ m " 1) (s,i,A)s k dsdA . 

T -co a(i) 

where 6(x) is the Dirac delta function. Performing the differentiation 
and using the relations H k (a(£),0] * 1, V[a(£),s,0] = 0, we reduce this 
equation to 

q k m) [aU),T] = 2u(£)6 kl 6(t) 

T 00 . (_) _ k - 

- f f 3 V [ a (£ ) , s ,x “A ] G, (s,£,X)s dsdA . (25) 

0 aU ) T k k 
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It follows from Eq. (23) Chat 

r k ' 1 x^ m) (r,£,T) = 2u(£)6 kl H k (r-,t) 

t X « ■ _ (■ 1 ) . k 

-f f. f 3 V [aU) ,s,A-A]G, " (s,£,A)H, (r,i-A)s dsdAdA 

0 0 .(» x k k k 

T 00 

+/ / V (r,s f T-A)G^ m ^ (s,2,A)s^dsdA . (26) 

0 a( £ ) * k 

Substituting this into Bq. (18), we have after some simple manipulation 

r k 1 (a) ( r ) _ y 2 u (£)3 H (r,£,T-£)d£ <5, , 

k t k Ki * 

-/ / / / 3 V [a(£),s,A-X]G^ ra " 1) (s,£,A)3 li ( r , £ , T - £- A MsdAdAd 

00 0 a(£) A k * 

T t-£ 00 / i \ k 

+/ / / 3 V, (r,s,T-£-A)G, '(s,£,A)s dsdAd£ . (27) 

0 0 aU) T k k 

With xl m ^(r,T) given by Eq. (27), the m-ch iterative solution given by Eqs . (6) 

1C 

and (13) is now complete. The convergence of the iterative .solution (6) 

k- 1 

may be proved by showing that there exisCs a time t below which (ek) 

c 

^ £ M, M being a bounded constant, for given U(t) and a(x). 


The extension of the present method to non— symmetric flow is immediate. 
In place of Eq. (6), the solution now can be written as 


m 


- -aa6 + Z (cR)^ sin k0 

k=l k 


m 


+ Z (sR) k 1 cos k0 

k-1 


’he Jacobian in Eqs. (7) to (11) will have to be modified to include cos k0- 


terms. A set of integral transforms similar to that given by Eq. (13) will 

then reduce the partial differential equations for and into a 

set of 2k -d imens Lona 1 heat equations. The solution of and bv 

the Laplace transform method will follow the line described for the solution 
c (m) 

or li’k given above.' 


* — 
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IV. Separation Point 


An appropriate definition of separation points in unsteady two- 

dimensional flows is of great current interest. For the discussions on 

the controversy and the relevant references, we refer the readers to the 
3 

work of Ho . Here, we are concerned with the separation points on a 
radially deforming cylindrical surface. With respect to this surface 
the velocity components (v r , ,v Q l ) are given by 

V = v r ' a • V = v 9 • 

If we define a separation point at such a surface to be the point where 
the flow relative to the moving surface just starts to reverse its direction, 
then we must have 


9 v a ' * 0 at r a a(t ) 
r 9 


Since v ' * v , this condition can be written as 


(9 <p) = 0. 

rr r=a 


It follows from Eq. (6) that at the separation points we must have 


m 


I ^(9 sin k9(i) a 0 , 

k-1 ” k r=a 


(29) 


where, with the aid of Eq. (13), 


.(n») x _ 2 u(t) „ k-2 (m)/_ , _ k "l/„ ..(m^ 

(d ip, ) = — r — 6, , - a x, (a,t) + a (9 X, ) _ 

rr k r=a a(t) kl k r k r=a 


(30) 


It is easily verified that at r=a(t), the last two integrals in Eq. (27) 
vanish, and Eq . (27) gives 


a k 1 x^ m ^ [a(t ) , t ] = 2u(t)6^ k 


as is required by the boundary condition Eq . (15). Thus the first two terms 

in Eq . (30) cancel each other. The last term in Eq . (30) can be obtained from 


L 3 


Eq . (27) by a simple differentiation. It can be shown that the shear 
stress also vanishes at the separation points for this particular problem. 

For the ease of comparisons with some known results, we define a 
new Reynolds number based on the uniform velocity Q = and the radius 
a o * i • e • > 

(u)6)a x toa 

R e ^ - ($-)(-?-) - eR. 

e v a v 

o 

Combination of Eqs. (29) and (30) then gives 

Z (R a) k ‘ 1 O r x, (m) ) sin ke = 0 . (31) 

k=l r=a 

Several comments concerning the general feature of two-dimensional flow 
separation can now be made. In the limiting case of creeping flows R e -*0, 
and Eq. (31) is reduced to 

( 3 r X^) sin 0 = 0. 

r=a 

Thus the only points on the cylinder where shear stress vanish are at 

the forward and rear "stagnation" points respectively at and 9=0. 

For a moderately large R such that the second term in Eq. (31) is comparable 

e 

to the first one, only the second iteration may be adequate to approximate 
the separated flow. Eq. (31) then gives 


.in 9 (3 x + 2R 3_x 0 cos 9) = 0 

*■ 1 e r 2 

r=a 


(32) 


Thus, in addition to the "stagnation" points there will be two symmetrically 

positioned primary separation points given by 

( 2 ) 


3 r X 1 

cos 9 = - ( rr s ■) 


2R e X 


( 2 ) 


(33) 


e ^ 2 


r=a 


bq . (33) has a solution only after a finite separation time is reached when 

the quotient in Eq . (33) becomes less than one in absolute value. At higher 
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Reynolds numbers, higher harmonic terms in 9 must be retained. This may 
allow secondary and further separation points to be obtained as the roots 
of the transcendental equation of order higher than two. It is obvious 
that separation may occur in mid-stream. This can be seen easily from 
Eq. (29) when r is evaluated at a point inside the flow rather than at r=a. 

It should be pointed out that if vortex shedding takes place, the flow 
is no longer symmetric, and separated flow must be based on Eq. (28). 
Asymmetric flow separation is not included in. this study. 

V . Results 

To* obtain numerical results, we used short-time expansions of Bessel 

functions appearing in the integrands of the analytical results given 

in the previous sections. For the evaluation of the multiple integrals, 

14 

we used the M-point Gauss quadrature formula . For 

each integration interval, the value of M was increased until the values 
of the given integral corresponding to the two successive values of 
M differ from each other by less than 100th of a percent. All numerical 
computation was carried out with double precision on an IBM 4341. 

To compare our iterative solution with some known experimental 
and theoretical results, we obtained numerical results for the case of 
a rigid cylinder first. The time dependent angle of separation was 
determined from Eq. (33) for the case of an impulsively started uniform motion 
of a rigid cylinder. The results are given in Fig. 2 together with 
other known results. The finite separation times for various values 
of Re are all slightly larger than those given by Collins & Dennis^ 
and Thoman & Szewezyk^. The difference may be due to the different 
initial conditions us.ed . The initial flow obtained by Collins & Dennis was 
based on a bounda rv- 1 aye r approximation. Thoman & Szewczyk and others who 
used pure numerical solutions assumed either that the flow was initially 
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irrotational everywhere or at some finite distance from the cylinder. The 

initial condition used in this study is that the fluid is completely quiescent 

everywhere. The separation angle for a given R g increases rapidly after the 

onset of separation and approaches an almost steady value within a small 

fraction of viscous diffusion time. Note that the dimensionless time t^ 

of the quoted theoretical studies are related to ours by x=t^/R e . Th e 

authors are not aware of any measurements of separation angles for the initial 

time smaller than x=0.1. We extrapolate the measurements of Coutanceau and 

Bouard^ for Re=20 in a tank whose diameter is 8.3 times that of the test 

cylinder. Our predicted angle is considerably larger than their measured value 

This is quite to be expected, since the wall of the test tank was shown by 

them to have the effect of reducing the separation angles. There is also 

considerable uncertainty involved in the extrapolation. It should be pointed 

out that the non-dimensional time t* of Coutanceau & Bouard is related' to 

ours by x=2t*/R . 

e 

Figure 3 demonstrates the effect of .the radial deformation of the 
cylinder, whicft~impulsively starts a constant velocity translation, on 
the separation angle. It is seen that, for a given R e , the radial expansion 
has an effect of reducing the onset time of separation and increasing the 
separation angle at the same instant. The larger the rate of expansion 
the larger is the effect. The earlier separation and larger separation 
angle are all associated with larger deceleration of fluid induced by 
radial expansion. Figure 4 shows Chat radial contraction has an 
opposite effect. Figure 5 shows the same effect demonstrated in Figs. 3 
and 4 except that the cylinder translation is now not an impulsively started 
uniform motion but is given by a constant acceleration followed by a 
constant velocity after some finite time x . This latter cylinder 
motion is obviously more easily attainable in experiments. 
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VI. Discussion 


An equation which describes the evolution of the separation angle on 
a deforming cylinder in an unsteady cross- flow was derived from an 
iterative solution of the Navier-Stokes equation for finite Reynolds 

numbers. The validity of the iterative solution hinges on the condition 
that the penalty functions remain negligible. This can be verified 
aposteriori. These penalty functions were found to be negligible in a 
larger period of time for the case of constant acceleration than for the 
case of impulsively started motion. Thus the time periods during which 
the iterative solution is valid can be ascertained aposteriori at each 
iteration. The advantage of the present method over other known methods 
is that only one computer program needs to be written for any u(x) and a(r). 
The solutions of the evolutin equation of separation point on a 
circular cylinder reveal that a radial contraction retards separation 
due to the enhanced acceleration along the wall. The converse is true 
for the case of radial expansion. In the limit of R^-K), It was shown 
that no separation at the wall other than that at the rear "stagnation" 
point can take place. All these results suggest that the flow separation 
on the body of revolution translating at an angle of attack in a fluid does 
nto originate at the nose-tip of the body, but originates somewhere downstream 
from the tip of the body. The results reported in this work were based on the 

second interative solution. Hence only the primary flow separation points 

were found . For the follow-up work, we intend to investigate the onset and 
evolution of the symmetric secondary flow separation at the deforming wall 
as well as in the mid-stream and the growth of the closed-wake, by use of the 
third iterative solution. The iterative solution of the form given by Eq. 

(28) for asymmetric flows will also be developed to study the breakup of 
a closed-wake and the ensuing vortex shedding. 
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Appendix 


A General Method of Superposition for Problems Involving 

a Moving Boundary 


We wish to solve for the unknown function F from the following 
partial differential equation 


( 3 t -L)F( £ ,t ) = G(£,t) , (A 


with the boundary 

condition 



2u(t) , 

T > 0 

F [a(i ) , t ] = 

{ 



0 

T < 0 


and the initial condition 

F(r,t) = 0 , t<0 , (A 

where t is time, r is a position vector, G(r,x) is a given function, 

L is any linear differential operator in space variables, and r=a(T) 
specifies the instantaneous position of the moving boundary at 
any time t . 

It is presently shown that the solution of this system is 
given by 


T 

- F ( r,"r )-=—/- — -3 " Ft? , T-Ttar - CX 

o T 

where F(r,£,i) is the solution of (A-l) subject to a fixed boundary _ 
condition 

1, t > 0 

F [ a( £, ) , 2.", t ] = 2u(DH(t) , H(r) = { , (A 

0, t < 0 


and the inicial condition 


i) 


2 ) 


3) 


■51 


L8 


F(r,£,T ] = 0 , 


T < 0 


( A-6 ) 


where l is a constant parameter. Since the differential system (A-l), 

(A-5) and (A-6) is invariant with respect to time translation, 

F(r,!,T-l) and F(_r,l,T-l-dl) 

are also solution of (A-l). Each satisfies (A-5) and (A-6) with t 
in both of them replaced respectively by t-1 and t-l-dl. It follows 
that dF given by 

dF * F(r,l,T-l) - F(£, 1 , t - 1-dl) 

» 3 F(r,l,T-l)dl (A-7) 

t — 

satisfies (A-l), since it is merely a linear superposition of two 
solutions of (A-l). Moreover (A-7) satisfies the boundary condition 

dF[a(l) ,1 ,t- 1] - 2u(l)[H(T-l)-H(T-l-dl)] (A-8) 

and the initial condition 

dF = 0 , t<1 . ( A-9 ) 

The integral of dF from 1=0 to 1 =t , i.e. 

T 

FtryrT — 3 F('r CA-lOy 

O 

being merely a linear superposition of solutions of (A-l), each satisfying 
(A-8) instantaneously from £,=0 to l-x for any x>0, is itself a solution 
of (A-l). Moreover it satisfies the boundary condition (A-2), since the 
integration of (A-8) gives 

T T 

/ dF[a(0,9..T-Md?.«Ffa(0,x] = / 2u ( l ) [ H( r - l )-H(r - l-d Z ) ] d tf = 2u(: ) 
o 0 

It also follows from (A-9), (A-10) and the fact that Z >0 , the initiaL 
condition ( A - 3 ) is satisfied by (A-10). 
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Note that (A-2) is the Dirichlet boundary condition. It is obvious 
that the same method applies to the Neumann or mixed boundary-value 
problems when the boundary conditions (A-2) and (A-5) are replaced by 
byoundary conditions specifying the spatial derivative of F or the 
linear combination of F and its spatial derivative. 
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Figure Captions 


Figure 1 
Figure 2 

Figure 3 

Figure 4 
Figure 5 


Definition sketch. 

Migration of separation point for the case of impulsively 
started uniform velocity. 

Effect of radial expansion on separation angle for the case of 

impulsively started uniform velocity. Rigid cylinder, . 

Deforming cylinders , a(x )=i+0.01x ; , 

a(r )=l+0.02x . 

Effect of radial contraction on separation angle. Rigid cylinder 

. Deforming cylinder: , a(x )=l-0. lx ; , 

a(x )®l-0.02x . 

Effects of acceleration of a deforming cylinder on separation 
angle; a(x ) a l+0. Olx . — -- , x 3 0.1; , x a 0.01. 
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